function moments_NLI = NLI_moment_funct(data,col_long)

% Calculate a vector of targeted moments related to the non-labor income process

NLI_orig = data(:,col_long.nonlabor_income);
NLI = NLI_orig;
dad_age = data(:,col_long.dad_age);
dad_age(dad_age<25 | dad_age>60) = NaN; % discard all men younger than 25 or older than 60 (too few observations)
dad_educ = data(:,col_long.dad_educ);
dad_educ(dad_educ<12) = NaN; % discard all men with less than HS degree (too few observations)

index_discard = (dad_age>0 & dad_educ>0);
NLI(index_discard == 0) = NaN; % discard these NLI observations

dummy_pos = NaN(size(data,1),1);
tmp0 = find(NLI==0);
tmp1 = find(NLI>0);
dummy_pos(tmp0) = 0;
dummy_pos(tmp1) = 1;
dummy_pos(index_discard == 0) = NaN;
% Note: this dummy is a NaN if NLI = NaN.

NLI_pos = NLI;
NLI_pos(NLI == 0) = NaN;
NLI_pos(index_discard == 0) = NaN;
% this vector only contains the strictly positive values of NLI

% List of moments we want to capture
moments_NLI = NaN*ones(21,1);

moments_NLI(1) = nanmean(NLI_pos);
moments_NLI(2) = nanstd(NLI_pos);
moments_NLI(3) = nanmean(NLI_pos .* dad_age);
moments_NLI(4) = nanmean(NLI_pos .* dad_age .* dad_age);
moments_NLI(5) = nanmean(NLI_pos .* dad_educ);
moments_NLI(6) = nanmean(NLI_pos .* dad_educ .* dad_educ);
moments_NLI(7) = nanmean(NLI_pos .* dad_age .* dad_educ);

% add some moments to capture the proportion of zeros better
moments_NLI(8) = nanmean(dummy_pos); % proportion of zeros
moments_NLI(9) = nanstd(dummy_pos);
moments_NLI(10) = nanmean(dummy_pos .* dad_age);
moments_NLI(11) = nanmean(dummy_pos .* dad_age .* dad_age);
moments_NLI(12) = nanmean(dummy_pos .* dad_educ);
moments_NLI(13) = nanmean(dummy_pos .* dad_educ .* dad_educ);
moments_NLI(14) = nanmean(dummy_pos .* dad_age .* dad_educ);

moments_NLI(15) = nanmean(NLI);
moments_NLI(16) = nanstd(NLI);
moments_NLI(17) = nanmean(NLI .* dad_age);
moments_NLI(18) = nanmean(NLI .* dad_age .* dad_age);
moments_NLI(19) = nanmean(NLI .* dad_educ);
moments_NLI(20) = nanmean(NLI .* dad_educ .* dad_educ);
moments_NLI(21) = nanmean(NLI .* dad_age .* dad_educ);

beg_age_dad = [25 34 39 44];
end_age_dad = [33 38 43 60];

tmp = NaN(4,3);
for t = 1:4
    % select the rows we need for this age group
    row_ind_dad = (data(:,col_long.dad_age) >= beg_age_dad(t) & data(:,col_long.dad_age) <= end_age_dad(t));
    tmp(t,1) = nanmean(NLI(row_ind_dad));
    tmp(t,2) = nanmean(NLI_pos(row_ind_dad));
    tmp(t,3) = nanmean(dummy_pos(row_ind_dad));
end
tmp = tmp(:); % vectorize
moments_NLI = [moments_NLI; tmp];

mom_age = data(:,col_long.mom_age);
mom_age(mom_age<25 | mom_age>55) = NaN; % discard all women younger than 25 or older than 55 (too few observations)
mom_educ = data(:,col_long.mom_educ);
mom_educ(mom_educ<12) = NaN; % discard all women with less than HS degree (too few observations)

tmp = NaN(12,1);

tmp(1) = nanmean(NLI_pos .* mom_age);
tmp(2) = nanmean(NLI_pos .* mom_age .* mom_age);
tmp(3) = nanmean(NLI_pos .* mom_educ);
tmp(4) = nanmean(NLI_pos .* mom_age .* mom_educ);

tmp(5) = nanmean(dummy_pos .* mom_age);
tmp(6) = nanmean(dummy_pos .* mom_age .* mom_age);
tmp(7) = nanmean(dummy_pos .* mom_educ);
tmp(8) = nanmean(dummy_pos .* mom_age .* mom_educ);

tmp(9) = nanmean(NLI .* mom_age);
tmp(10) = nanmean(NLI .* mom_age .* mom_age);
tmp(11) = nanmean(NLI .* mom_educ);
tmp(12) = nanmean(NLI .* mom_age .* mom_educ);

moments_NLI = [moments_NLI; tmp];

beg_age_mom = [25 33 38 43];
end_age_mom = [32 37 42 55];

tmp = NaN(4,3);
for t = 1:4
    % select the rows we need for this age group
    row_ind_mom = (data(:,col_long.mom_age) >= beg_age_mom(t) & data(:,col_long.mom_age) <= end_age_mom(t));
    tmp(t,1) = nanmean(NLI(row_ind_mom));
    tmp(t,2) = nanmean(NLI_pos(row_ind_mom));
    tmp(t,3) = nanmean(dummy_pos(row_ind_mom));
end
tmp = tmp(:); % vectorize
moments_NLI = [moments_NLI; tmp];

end